Fractional charge perspective on the band-gap in density-functional theory 
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The calculation of the band-gap by density-functional theory (DFT) methods is examined by 
considering the behavior of the energy as a function of number of electrons. It is found that the 
incorrect band-gap prediction with most approximate functionals originates mainly from errors in 
describing systems with fractional charges. Formulas for the energy derivatives with respect to 
number of electrons are derived which clarify the role of optimized effective potentials in prediction of 
the band-gap. Calculations with a recent functional that has much improved behavior for fractional 
charges give a good prediction of the energy gap and also £homo — — / for finite systems. Our results 
indicate it is possible, within DFT, to have a functional whose eigenvalues or derivatives accurately 
predict the band-gap. 
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One of the many important uses of density-functional 
theory (DFT) is the calculation of the band-structure 
which has many applications throughout physics, for ex- 
ample in semiconductors, electron transport and reac- 
tions at surfaces. The first step in achieving accuracy in 
the band structure is to understand the band-gap which 
standard functionals have long been known to systemati- 
cally underestimate by as much as ^50%. Recent efforts 
have focused on use of the optimized effective potential 
(OEP) method, which can often give an improvement in 
the prediction of band-gaps for small-gap semiconduc- 
tors, but has problems with wider gap semiconductors 
and insulators [TJ [3J. In many cases it has proved 
necessary to move to the quasi-particle GW theory to 
calculate the band-gap of solids accurately [4] . Conven- 
tionally, the band-gap problem, has been related to the 
so-called "derivative discontinuity" in the exchange corre- 
lation potential: even with an accurate Kohn-Sham po- 
tential, the energy-gap is still different from the true gap 
by an amount of the derivative discontinuity [5 , 6 . This 
perspective, however, does not offer the understanding or 
the mechanism needed for making progress for band-gap 
prediction with DFT. 

In this Letter, a new perspective is offered: the band- 
gap problem is shown to be related to the behavior of 
approximate density functionals for fractional numbers 
of electrons, an issue which has drawn considerable re- 
cent interest [71 [51 IS]- This enables us to understand 
the problem with band-gap calculations and offers ideas 
to develop functionals which predict the band-gap cor- 
rectly. Examples will be given for molecules where the 
energy gap can be compared with explicit calculations of 
systems with fractional charges. 

The fundamental band-gap for an iV-electron system 
in an external potential v(r), is given by 

E g = [E v (N-l)-E v (N))-[E v (N)-E v {N + l)) 
= I -A, (1) 

where E V (N) is the ground-state energy of the N par- 



ticle system and I is the ionization energy and A is the 
electron affinity. For a system with a fractional num- 
ber of electrons N + SN, with < SN < 1, it has been 
shown that the energy is a straight line connecting the 
total energies at integer numbers of electrons; namely, 
E V (N + SN) = SNE V (N) + (1 - SN)E V (N +1) PHE]. 
This linear relation means that the energy gap in Eq. 
can be given by the derivative difference 



lim • 

SN^O 



dE v 
dN 



dE„ 



N+SN 



dN 



(2) 



N-SN 



If we substitute in the DFT total energy expression E v 
T s [p] + Vextlp] + J[p] + E xc [p], we have 



E„ 



(3) 



where the E^ s is the gap in a Kohn-Sham calculation 
and the A xc is the derivative discontinuity [SJ[5]. 

In this work we identify the problem with calculations 
using approximate density functionals by considering the 
basic assumption in Eq. ([2]), that the energy at N + 1 
and at N — 1 can be given simply from the derivatives 
at N. This is true for exact DFT but it may or may 
not be true for approximate functionals. The key is to 
investigate the behavior of the total energy as a function 
of numbers of electrons. To do this we consider a non- 
interacting ground state reference system where we allow 
the occupation numbers of the orbitals to vary the num- 
ber of electrons smoothly; the first-order reduced density 
matrix of the reference system is given by 



P s (r,r') 



5><Mr)0*(r') 



(4) 



where rii = 1 for i < f, rii = SN for i = f, and rii = 0, 
for i > f, and f is the index for the frontier orbital. 
The electron density is just its diagonal p s (r) = p s (r,r). 
The behavior of three qualitatively different exchange- 
correlation functionals: the local density approximation 
(LDA), Hartree-Fock (HF) and MCY3 [H] is shown in 
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Figure 1: Behavior of the energy of the carbon atom with 
number of electron with fractional charges, for several differ- 
ent functionals minimized with and without OEP. Dotted line 
follows the initial slope for the non-straight functionals. Inset 
shows 6 < TV < 7 range in more detail. 

Fig. 1. MCY3 was constructed as a functional of p s (r, r') 
to give a much improved description of fractional num- 
bers of electrons and we can see that it gives a straight 
line interpolation between the integers. The interpola- 
tion from the other functionals is far from linear, HF 
curves in a concave manner and LDA in a very convex 
manner. There are two main problems with approximate 
density functionals: First, they can have a remarkably 
different behavior from the exact functional in fractional 
charge systems, compared to their behavior for the corre- 
sponding integer charge systems. Second, the error in the 
integer charge systems can also be significant. The com- 
bined effects lead to the error in the band-gap prediction 
from derivative information. 

For molecules, LDA has a very reasonable description 
of the integer values (/ and A are given well) but a much 
worse description in between the integers. The use of the 
first derivatives for LDA will clearly not give the / and A 
from the integer calculations. Because of the convexity 
of the curve, LDA will give too small a value for / and 
too large a value for A, meaning that the band-gap I — A 
will be too small as shown by the dotted lines in Fig 1. 
The case for HF is very different as the integer values are 
not as good because of the lack of correlation, and also 
it curves in a concave manner. These two errors cancel 
each other to some extent in the prediction of / but add 
together in the calculation of A. For HF the value of I 
will be about right and A too small, meaning that the 
band-gap I — A will be too large as shown in Fig 1 . For 
MCY3, as it is very straight, the use of the derivatives 
will give a prediction very similar to the integer calcu- 
lations for / and A, and the band-gap I — A should be 
accurate as the integer values are good. For functionals 
that have a linear behavior for fractional charge systems 
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Figure 2: Comparison of MCY3 eigenvalues from GKS and 
OEP calculations and also including A xc . 

(e.g. MCY3 and the exact functional) the initial deriva- 
tive is all that is needed to calculate I and A and the 
band-gap as in Eq. ([2|. Hence, we consider analytic 
expressions for dE v /dN for some families of exchange- 
correlation functionals. The main ideas and results are 
summarized here, with further details presented in the 
supplementary material [13 . 

In the fractional charge non-interacting system, Eq. 
Q, the orbitals are the eigenstates of an one-electron 
local potential v s (r) 

(-l^ + vs^j \<t>i) =£i|&), (5) 
or a non-local potential v^ L (r,r'). 

(-iv a +^)|^)=ef«|^>. (6) 

The former is the original Kohn-Sham (KS) reference 
system and the latter has been called the Hartree-Fock- 
Kohn-Sham (HFKS) [TJ or the generalized Kohn-Sham 
(GKS) method [15]. 

For the Kohn-Sham reference system with local po- 
tential v s (r), we here use the potential- functional for- 
mulation [IB]. The electron density p s (r) can be rep- 
resented as the set of orbitals and occupation numbers 
{<j>i,rii}, or equivalently as the local KS potential and 
total particle number {v s (r), N}. Thus the total en- 
ergy functional, formally in term of density as E v [p s (r)], 
can be expressed as E v [v s (r) , N] . The ground state en- 
ergy is the minimum of the KS energy functional, ex- 
pressed (explicitly or implicitly) in terms of the local po- 
tential i> 8 (r): E V (N) = min^ E v [v s , N] = E v [vf s , N] , 
where the minimizer vf s is the optimized effective po- 
tential (OEP), as established recently [IBJ. The vari- 
ational nature of uf s simplifies the calculation of the 
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Table I: Comparison of et against corresponding experimental numbers for LDA, HF and MCY3 for more details see [13] 



flumo Shomo I A 

Mol LDA HF MCY3 Expt 


Ehomo 

LDA HF MCY3 Expt 


LDA HF MCY3 Expt 


C 0.08 12.76 10.03 10.00 
O 0.23 16.80 11.57 12.16 
F 2 3.38 20.44 13.75 14.40 
OH 0.10 16.56 11.23 11.40 
MAE 4.678 2.082 0.319 


6.09 11.94 11.12 11.27 
7.28 14.11 13.01 13.62 
9.53 18.13 15.17 15.70 
7.21 13.90 12.70 13.20 
5.059 0.990 0.484 


6.01 -0.82 1.09 1.27 
7.05 -2.70 1.44 1.46 
6.15 -2.31 1.42 1.30 
7.11 -2.67 1.47 1.80 
4.295 3.328 0.488 



Errors from the explicit calculation of / and A |13| 



MAE 0.920 1.248 0.789 0.659 1.235 0.505 0.328 1.655 0.475 



j ■ ,. dEv(N) dEv[v g °,N] 

derivative: ^ = ^ — 



. Consider a change 

in the total number of electrons N = Nq + SN, where 
A*o is an integer and |5JV| < 1. At fixed vf*, all the or- 

bitals j</^ s |, as its eigenstates, are fixed. Since p s (r) 

is the ground state density of the reference potential vf s , 
only the frontier level occupation rif is allowed to change 
SN = 5n { , thus 



dE v (N) _ ( dE v [{<fi' ,m}] 



ON 



dn f 



(7) 



where the frontier orbital is either the LUMO, nf — 
niumojif SN > , or the HOMO, n f = ra ho mo,if SN < 0. 
We consider three cases for which the analytic derivatives 
can be obtained [T3"] : 

Case A: E xc [p s (r)], an explicit functional of p s (e.g. 
LDA or GGA): 



dE v (N) 
dN 



(8) 



where Sf is the KS eigenvalue for the frontier orbital in 
the local potential v s (r) — v(r) + v,j(r) + v xc (r). This is 
exactly the combination of Eq. ([7]) with Janak's theorem 
for m [17]. 

Case B: E xc [p s (r, r')], a functional of the first order 
density matrix minimized with a local potential as in Eq. 
^ (e.g. OEP exact exchange, EXX). 



dE v (N) 
dN 



= £f + (<j)f\v + VJ + V. 



NL 

xc 



V s \(j>f 



(9) 



where the non-local potential v xc L (r, r') = ^^%^y - ■ 
Eq. ^ is a key result, showing that for general orbital 
functionals, is not given by the frontier OEP eigen- 
value, £f , but with a correction term A xc . This general 
result agrees with [T] in the case of exact exchange, and 
is related to the results of [TB] from the self energy. 

Case C: E xc [p s (r, r')], with the energy minimized with 
respect to the orbitals <pi (e.g. HF ): 



9E V (N) _ GKS 

~dN~- e * 



(10) 



where ef K s is the eigenvalue of the frontier orbital of the 



non-local potential v + vj 



,NL 



(r,r'), as in Eq. m 



All three cases can be unified in the expression 



dE v {N) 
dN 



cll\ 



(11) 



where H Q g — — |V 2 



< + vj + v xc (r) for Case A where 
E xc = £Up(r)],"and H cS = ~\V 2 + v + v.j + v* c L {r, r') 
for case B and C where E xc = E xc [p s (r', r)], which are 
often called orbital functionals. In cases A and C, dE g^ 
is equal to the corresponding eigenvalue, but not in case 
B. The only difference in between cases B and C is the 
orbitals used to evaluate the overall expression. 

Combining Eqs. ^ and ( 11 ) thus leads to the general 
expression for the band-gap from derivatives for an N- 
electron system: 



E] 



((t>N+l\El c g\(j)N^ 



^n\H cB \0n)- (12) 



We consider a few illustrative atomic and molecular 
systems for which we have performed self-consistent cal- 
culations using a cc-pVQZ basis set in an modified ver- 
sion of CADPAC. We compare — £homo with the experi- 
mental /, — £iumo with the experimental A and also their 
corresponding differences. The results for LDA, HF and 
MCY3 are given in Table 1 . MCY3 gives very good agree- 
ment between — £homo and I which is to be expected from 
its straight line behavior. We should emphasize that this 
has not been seen before for calculations with approxi- 
mate exchange-correlation functionals. The error is rel- 
atively small, 0.5eV, and is roughly similar in I, A and 
also the difference I — A. LDA eigenvalues have a large 
error, with a consistent underestimation of I by about 5 
eV and overestimation of A by about 4 eV and a poor 
predicition of the gap. LDA does well for the explicit 
calculation of the N + 1 and N — 1 systems, it is just the 
use of the derivatives at N that lead to large errors, this 
is now clearly understood from the fractional charge pic- 
ture. The HF — £homo is often close to /, however there 
are larger errors for A and also for the gap. 

The band-gap issue is well understood for calculations 
with local density functionals (case A) or orbitals func- 



4 



Table II: Comparison of GKS and OEP eigenvalues, and 
dE/dn t for HF and MCY3, see p3]. 



for systems with a fractional charge. We have recently 
developed a functional with this linear behavior giving 



Mol 


MCY3 


MCY3 


ML Y o 


HF 


HF 


xlr 


Expt 




OEP 


vi K b 


OEP 








e 


dE/dni 




e 


dE/dnt 






C I-A 


0.70 


10.16 


10.03 


1.47 


13.49 


12.76 


10.00 


I 


10.58 


11.11 


11.12 


11.97 


11.94 


11.94 


11.27 


A 


9.88 


0.96 


1.09 


10.50 


-1.55 


-0.82 


1.27 


F 2 I-A 


4.06 


13.74 


13.75 


5.62 


20.49 


20.44 


14.40 


I 


14.67 


15.16 


15.17 


15.94 


18.11 


18.13 


15.70 


A 


10.61 


1.42 


1.42 


10.32 


-2.37 


-2.31 


1.30 



/ and — £i u 



A and a good prediction of 



tionals (case C). We now consider case B, orbital func- 
tionals in an OEP calculation, using the Yang-Wu direct 
minimization method [13120]. In Fig. 1 the OEP mini- 
mized energy is remarkably similar to the GKS minimized 
energy in both integer and fractional charge systems. We 
would therefore expect dE / 8Nn±$n to be the same as 
the GKS derivatives. 

Table 2 shows the eigenvalues from an OEP calcula- 
tion using a Fermi- Amaldi base potential, which has the 
correct asymptotic behavior. The asymptotic form of the 
potential has a large effect on the OEP eigenvalues, but 
not on the eigenvalue differences or the energy deriva- 
tives. The OEP(MCY3) — £homo — I (as is proven for the 
exact functional however the E^ s is much smaller 
than the exact gap. The inclusion of A xc gives a much 
better agreement between the derivatives and the GKS 
eigenvalues. This brings us on to the nature of this term; 
it is only the difference between a KS and GKS calcula- 
tion and is needed to correctly give the derivative at N. 
It does not, however, address the question of whether the 
functional used for the calculation has the correct straight 
line behavior for fractional numbers of electrons, which 
is the key question in the evaluation of the band-gap. 

Fig. 2 shows the behavior of the eigenvalues for car- 
bon with different numbers of electrons using MCY3. 
The GKS £homo for a fractional system is almost con- 
stant between integers due to the straight line behavior 
of MCY3. The OEP e f are markedly different to the GKS 
ef, but they become almost identical upon inclusion of 
A xc , which can be understood from comparing Eqs. (|9j) 



and ([T0f. The LUMO at N - SN is connected to the 
HOMO at N + 5N. This is clear from Eq. Q as the 
change in the number of electrons, N, is only through 
the frontier occupation numbers, rif, and the potential 
and, therefore, cigcnfunctions remain fixed, i.e. there is 
no mysterious discontinuity in the eigenvalues. 

In conclusion we have carried out analysis and calcula- 
tions on systems with fractional numbers of electrons to 
gain understanding of the band-gap problem in DFT. We 
show that the band-gap is only given by the eigenvalue 
difference if the functional has the correct linear behavior 



the band-gap in molecules. We have also considered OEP 
calculations, in which the derivative of the energy with 
respect to number of electrons is not given by the OEP 
eigenvalue. When the derivative is correctly evaluated, it 
gives practically the same as in GKS calculations. Our 
work thus provides the new insight: it is possible to have 
a functional which gives the correct band-gap from the 
eigenvalues or derivative information, so long as it has the 
correct fractional charge behavior and accurate energies 
for integer systems. Such possible functionals include ex- 
plicit functionals of the electron density E xc [p(v)\. We 
have only considered the explicit calculation of molecules 
but the same ideas are undoubtedly of key importance 
in solids. The understanding gained in this Letter offers 
a new perspective and way forward for accurate calcula- 
tions of the band-gap in DFT. 

Support from NSF is greatly appreciated. 
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